.ipynb .pdf

Binder

Contents

Wannierization

Electronic simulation codes like Quantum Espresso calculate the Bloch bands \(\varepsilon_{\mathbf kn}\) and wavefunctions \(\psi_{n\mathbf k}(\mathbf r)\) of a crystal lattice. These bands are gauge invariant while the wavefunctions allow for an ambiguity in the choice of complex phase. These gauges need to be consistent if we wish to transform these wavefunctions to a localized Wannier function :cite:`wannier_review`:

(3)\[ |\mathbf Rn\rangle=\frac{V}{(2\pi)^2}\int_{BZ}\text d\mathbf k\, \text e^{-i\mathbf k\cdot\mathbf R}|\psi_{n\mathbf k}\rangle. \]

Projection method

One way to remove the gauge inconsistency is to project trial orbitals \(|g_n\rangle\) on the basis of the Bloch wavefunctions:

(4)\[ |\phi_{n\mathbf k}\rangle=\underbrace{\sum_m|\psi_{m\mathbf k}\rangle\langle\psi_{m\mathbf k}}_{P_\mathbf k}|g_{n}\rangle. \]

The projector \(P_\mathbf k\) is gauge invariant so the \(|\phi_{n\mathbf k}\rangle\) are as well. The \(|\phi_{n\mathbf k}\rangle\) are however not necessarily orthogonal. In general we have:

(5)\[ \mathbf S_{mn,\mathbf k} = \langle\phi_{m\mathbf k}|\phi_{n\mathbf k}\rangle\neq\delta_{mn} \]

Marzari et al. :cite:`wannier_review` propose to orthogonalize the states using a Löwdin transform which finds a set of orthonormal basis functions \(|\tilde\psi_{n\mathbf k}\rangle\) that minimizes:

(6)\[ \sum_n \int\text d \mathbf x\,\left |\phi_{n\mathbf k}(\mathbf x)-\tilde\psi_{n\mathbf k}(\mathbf x)\right|^2, \]

which corresponds to taking:

(7)\[ |\tilde\psi_{n\mathbf k}\rangle=\sum_m|\phi_{m\mathbf k}\rangle \left(\mathbf S_\mathbf k^{-1/2}\right)_{mn}. \]

Example: 1 band

We illustrate the projection method using an example of a 1d isolated band with mock Bloch wavefunction \(\psi_\mathbf k(\mathbf r)=u_\mathbf k(\mathbf r)\text e^{i\mathbf k\cdot\mathbf r}\). This functions must obey the periodicity requirements:

• \(\psi_{\mathbf k+\mathbf G}(\mathbf r)=\psi_{\mathbf k}(\mathbf r)\), with \(\mathbf G\) a reciprocal lattice vector \(\rightarrow u_{\mathbf k+\mathbf G}(\mathbf r)=\text e^{-i\mathbf G\cdot\mathbf r}u_{\mathbf k}(\mathbf r)\).

• \(u_\mathbf k(\mathbf r+\mathbf T)=u_\mathbf k(\mathbf r)\), with \(\mathbf T\) a lattice vector.

Since \(u_\mathbf k(\mathbf r)\) is periodic in real space we can express it as a Fourier series:

(8)\[ u_\mathbf k(\mathbf r)=\sum_\mathbf G C_{\mathbf G+\mathbf k}\text e^{i\mathbf G\cdot r} \]

with the \(C_{\mathbf G+\mathbf k}\) some Fourier components and the sum over all reciprocal lattice vectors. We show that this choice complies with the periodicity of the Bloch function in \(\mathbf k\) space:

(9)\[ u_{\mathbf k+\mathbf G'}(\mathbf r)=\sum_\mathbf G C_{\mathbf G+\mathbf G'+\mathbf k}\text e^{i\mathbf G\cdot r} =\sum_\mathbf {\tilde G}C_{\mathbf {\tilde G}+\mathbf k}\text e^{i(\mathbf {\tilde G}-\mathbf G')\cdot r}=\text e^{-i\mathbf G'\cdot r}\sum_\mathbf {\tilde G}C_{\mathbf {\tilde G}+\mathbf k}\text e^{i\mathbf {\tilde G}\cdot r}=\text e^{-i\mathbf G'\cdot r}u_{\mathbf k}(\mathbf r) \]

Now we plug this form of \(u\) in equation numref{eq:wannier_transform} to get the Wannier function for \(R=0\):

(10)\[ w_0(\mathbf r)\propto\int_{BZ}\text d\mathbf k\, \text e^{i\mathbf k\cdot\mathbf r}\sum_\mathbf G C_{\mathbf G+\mathbf k}\text e^{i\mathbf G\cdot r} =\sum_\mathbf G\int_{BZ}\text d\mathbf k\, C_{\mathbf G+\mathbf k}\text e^{i(\mathbf G+\mathbf k)\cdot r}=\int\text d\mathbf k\, C_{\mathbf k}\text e^{i\mathbf k\cdot r}, \]

where the last integral is over all \(\mathbf k\) space.

This description makes clear that \(C_k\) needs to have power law or Gaussian decay to acquire exponentially localized Wannier functions. If we neglect this condition the according Bloch waves become highly irregular and oscillate wildly in the unit cell.

For this example we choose:

\[ C_\mathbf k = \text e^{-k^2/100}\cos(k/3.4), \]

which leads to the following Bloch function:

We note that the real part of \(\psi_k(x)\) is symmetric while the imaginary part is antisymmetric. Any symmetric \(g\) will then have a purely real overlap wth \(\psi\). For this function we choose a Gaussian:

(11)\[ g(x)=\text e^{-x^2}. \]

Any normalization of \(g\) and \(\psi\) will cancel in the construction of the orthonormal \(\tilde\psi\). We now do the inverse transform of the Bloch wavefunctions of ref{eq:wannier_transform} do get the following localized Wannier function:

Fig. 12 Three periodic \(u_i(\mathbf r)\) in a triangular lattice based on sp\(^2\) orbitals. We aim to construct localized Wannier orbitals using the three standard atomic orbitals: s, p\(_x\) and p\(_y\)

We see that the Wannier function is exponentially localized around \(x=0\) with local maxima (minima) appearing around \(x\in \mathbb Z\) which become exponentially smaller for larger \(x\).

Example: 2d atomic orbitals as Bloch wavefunction

To make matters a bit more interesting we’ll look at the hibridization of the s, p_\(x\) and p\(_y\) orbitals in a triangular lattice. As \(u\) component of the Bloch waves we take the three sp\(^2\) orbitals which we label with \(\gamma_m\) and as \(g_n\) we take the normal orbitals. We take the \(C_{i\mathbf k}\) to be the inverse Fourier transforms of the \(\gamma_m\):

(12)\[ C_{m\mathbf k}\propto\int\textbf d\mathbf r\,\gamma_m(\mathbf r)\text e^{-i\textbf r\cdot k}, \]

which we then plug in equation numref{eq:u_transform}:

(13)\[\begin{split} \begin{align} u_{m\mathbf k}(\mathbf r)&\propto\sum_\mathbf G\text e^{i\mathbf G\cdot r}\int\textbf d\mathbf {r}'\,\gamma_m(\mathbf{r}')\text e^{-i\mathbf{r}'\cdot (\mathbf G+\mathbf k)}\\ &=\int\textbf d\mathbf {r}'\,\gamma_m(\mathbf{r}')\text e^{-i\mathbf{r}'\cdot\mathbf k}\sum_\mathbf G\text e^{i\mathbf G\cdot (\mathbf r - \mathbf{r}')}\\ &=\int\textbf d\mathbf {r}'\,\gamma_m(\mathbf{r}')\text e^{-i\mathbf{r}'\cdot\mathbf k}\sum_\mathbf R\delta(\mathbf R+\mathbf r - \mathbf{r}')\\ &=\sum_\mathbf R\gamma_m(\mathbf{r}+\mathbf R)\text e^{-i(\mathbf{r}+\mathbf R)\cdot\mathbf k}, \end{align} \end{split}\]

or:

(14)\[ \psi_{m\mathbf k}(\mathbf r)=\frac1\Omega\sum_\mathbf R\gamma_m(\mathbf{r}+\mathbf R)\text e^{-i\mathbf R\cdot\mathbf k}. \]

with \(\Omega\) a normalization constant, where we follow the normalization convention of cite{wannier_review}:

(15)\[\begin{split} \begin{align} \frac{4\pi^2}{V}=\langle \psi_{m\mathbf k}|\psi_{m\mathbf k}\rangle &= \frac1{\Omega^2}\frac1N\sum_{\mathbf R\mathbf R'}\text e^{i(\mathbf R-\mathbf R')\cdot\mathbf k}\int\text d\mathbf r\,\gamma_m(\mathbf{r}+\mathbf R)\cdot\gamma_m(\mathbf{r}+\mathbf R')\\ &=\frac1{\Omega^2}\sum_{\mathbf R}\text e^{i\mathbf R\cdot\mathbf k}\int\text d\mathbf r\,\gamma_m(\mathbf{r}+\mathbf R)\cdot\gamma_m(\mathbf{r})\\ &=\frac1{\Omega^2}\int\text d\mathbf r\,\gamma_m(\mathbf{r})^2\\ &=\frac1{\Omega^2} \end{align} \end{split}\]

where the factor N is the amount of unit cells in the lattice which cancels with one of the summations. Summations over nonzero \(\mathbf R\) cancel due to lattice symmetry.

Fig. 13 Three periodic \(u_i(\mathbf r)\) in a triangular lattice based on sp\(^2\) orbitals. We aim to construct localized Wannier orbitals using the three standard atomic orbitals: s, p\(_x\) and p\(_y\)

We write ref{eq:nonorthogonal} out in matrix form:

\[\begin{split} \begin{pmatrix} |\phi_{1\mathbf k}\rangle\\ |\phi_{2\mathbf k}\rangle\\ |\phi_{3\mathbf k}\rangle\\ \end{pmatrix}=\underbrace{ \begin{pmatrix} \langle \psi_{1\mathbf k}|g_{1}\rangle & \langle \psi_{2\mathbf k}|g_{1}\rangle &\langle \psi_{3\mathbf k}|g_{1}\rangle\\ \langle \psi_{1\mathbf k}|g_{2}\rangle & \langle \psi_{2\mathbf k}|g_{2}\rangle &\langle \psi_{3\mathbf k}|g_{2}\rangle\\ \langle \psi_{1\mathbf k}|g_{3}\rangle & \langle \psi_{2\mathbf k}|g_{3}\rangle &\langle \psi_{3\mathbf k}|g_{3}\rangle\\ \end{pmatrix}}_{A_\mathbf k} \begin{pmatrix} | \psi_{1\mathbf k}\rangle\\ | \psi_{2\mathbf k}\rangle\\ | \psi_{3\mathbf k}\rangle\\ \end{pmatrix}, \end{split}\]

then \(S_\mathbf k=A_\mathbf k^\dagger A_\mathbf k\), which we diagonalize: \(S_\mathbf k=U_\mathbf k^\dagger S^d_\mathbf k U_\mathbf k\). Now we get: \(S_\mathbf k^{-1/2}=U_\mathbf k \left(S^d_\mathbf k\right)^{-1/2} U_\mathbf k^\dagger\) and \(S_\mathbf k^{-1/2}A_\mathbf k\) is the rotation matrix for the Bloch functions. We remind that any normalization of \(g_i\) cancels in this product.

The components of \(A_\mathbf k \) can be simplified analytically:

\[\begin{split} \begin{align} \langle \psi_{m\mathbf k}|g_n\rangle&=\int\text d\mathbf r\sum_\mathbf R\gamma_m(\mathbf{r}+\mathbf R)\text e^{i\mathbf R\cdot\mathbf k}g_n(\mathbf r)\\ &=\sum_\mathbf R\text e^{i\mathbf R\cdot\mathbf k}\int\text d\mathbf r\,\gamma_m(\mathbf{r}+\mathbf R)g_n(\mathbf r), \end{align} \end{split}\]

where the integral does no longer depend on \(\mathbf k\). We compute this integral using sums over coordinate meshes to speed up calculations.

Fig. 14 The three Wannierfunctions which resemble the \(g_i\): \(\tilde\psi_1\) takes after \(p_x\), \(\tilde\psi_2\) after \(p_y\) and \(\tilde\psi_3\) after \(s\).

The three Wannier wavefunctions in The three Wannierfunctions which resemble the g_i: \tilde\psi_1 takes after p_x, \tilde\psi_2 after p_y and \tilde\psi_3 after s. resemble the atomic orbitals chosen for the \(g_i\). The \(p\) orbitals delocalize to neighboring sites.

Wannierization Graphene